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Abstract 

Monte Carlo simulation has been performed in a two-dimensional modified XY- 
model first proposed by Domany et.al [E. Domany, M. Schick and R. H. Swendsen, 
Phys. Rev. Lett. 52, 1535 (1984)]. The cluster algorithm of Wolff has been used 
and multiple histogram reweighting is performed. The first order scaling behavior of 
the quantities like specific heat, order parameter susceptibility and free energy barrier 
are found to be obeyed accurately. While the lowest order correlation function was 
found to decay to zero at long distance just above the transition, the next higher order 
correlation function shows a non-zero plateau. 
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1 Introduction 

More than two decades ago Domany et.al [lj proposed a generalization of the two-dimensional 
XY- model where the shape of the usual cos# type potential could be modified with the 
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help of a single parameter. The two-dimensional spins located at the sites of a square lattice 
interact with the nearest neighbors through a potential 



where is the angle between the spins and p 2 is a parameter used to alter the shape of the 
potential. For p 2 = 1 the potential reproduces the conventional XY-model while for larger 
values of p 2 the potential well becomes narrower. The conventional two-dimensional XY 
model does not possess any true long range order which is ruled out by the Mermin Wagner 
theorem. However a continuous quasi-long-range-order-disorder transition resulting from 
the unbinding of topological defects (2J 13] is known to occur in this system and the order 
parameter correlation function is characterized by a slow algebraic decay instead of the fast 
exponential decay observed in a disordered system and this is referred to as the Kosterlitz- 
Thouless (KT) transition in literature. Domany et.al [TJ performed Monte Carlo (MC) 
simulation and observed that as the potential well gets narrower with the increase in the 
parameter p 2 , the continuous transition gets converted into a first order phase transition and 
for p 2 = 50 the transition is very sharp as is manifested by a huge peak in the specific heat. 
This phenomenon is in apparent contradiction with the prediction of the renormalization 
group theory according to which systems in the same universal class (having same symmetry 
of the order parameter and same lattice dimensionality) should exhibit the same type of 
phase transition with identical values of critical exponents. 

The generalized XY-model of Eqn. (TTJ has been analyzed by a number of authors [U |S] 
using the renormalization approach of the Migdal-Kadanoff type. These investigators were 
of the opinion that the transition in the generalized XY-model appears to be first order in 
nature because the MC simulation of Domany et.al [T] and Himbergen [B] were carried out 
on relatively small lattices and for large system sizes the usual KT transition is expected to 
occur. Nearly a decade later, Mila [7j using the same sort of renormalization group analysis 
arrived at a similar conclusion. Lastly, using the same line of approach, Garel et.al [8] put 
forward a different type of interpretation of the above mentioned RG analysis and were of 
the view that the transition is indeed first order. 

Minnhagen [9j [TOj [11] has carried out a detailed study of the behavior of the phase 
transition exhibited by a 2-D Coulomb gas, which very well describes the characteristics of a 
2-D system consisting of vortex-antivortex pairs. It was demonstrated that the KT behavior 
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is obtainable in a 2-D Coulomb gas only at low particle densities. For higher particle densities 
the charge unbinding transition was shown to be first order. Also a new gas-liquid like critical 
point was found in the 2-D Coulomb gas — the first order line in the temperature-particle 
density plane ends at a critical point. The KT transition line, obtainable at lower densities 
was seen to join smoothly with the first order line at a temperature slightly lower than 
the critical point. Jonsson, Minnhagen and Nylen [12] performed MC simulation in a 2-D 
XY-model with a modified potential, which essentially is equivalent to that of Eqn.fll]) and 
established a new critical point. They determined the critical exponents for the system and 
interpreted the transition to be of the vortex unbinding type. 

van Enter and Shlosman [13] presented a rigorous proof that various SO(n)-invariant 
n-vector models which have a deep and narrow potential well, would exhibit a first order 
transition. The model represented by Eqn. (JTJ) is a member of this general class of systems. 
These authors based their proof on the so called method of reflection positivity, a technique 
borrowed from the field theory and used in statistical mechanics, van Enter and Shlosman 
argued that in spite of the order parameter in 2-D n-vector model being predicted to vanish 
by the Mermin- Wagner theorem, long range order prevails in the system via higher order 
correlation functions. More recently, S. Ota and S. B. Ota [H] have performed MC simu- 
lation of the modified XY- model using microcanonical ensemble and have identified a first 
order phase transition in the system. 

The present article describes MC simulation of the 2-D modified XY-model where com- 
putations have been performed on systems of reasonably large size and finite size scaling 
rules for first order phase transition have been tested on the results of the simulation. The 
motivation is to resolve the question on the nature of the phase transition in this model 
and the contradictions among the views put forward by different investigators for the last 
quarter of a century as has been summed up above. Our observation is that the transition 
is indeed first order for a large value of the parameter p 2 (we have used p 2 = 50) as all finite 
size scaling rules are nicely obeyed. We however have made no attempt to investigate the 
existence of the critical point in this model or to determine the critical exponents as has 
been done by Jonsson et.al [12] in relatively small systems. Among other observables we 
have computed the spin-spin angular correlation functions of different orders. We observe 
that while the lowest order correlation function decays to zero, the next higher order corre- 



lation function has a finite plateau which is in accordance with statement of van Enter and 
Sholsman |13j . 

Another interesting aspect of our work is the application of the Wolff cluster algorithm 
[T5] to simulate the model. It has been pointed out by a numbers of workers [TJ [TB] that 
the two-dimensional model is difficult to simulate using the conventional single spin flip 
Metropolis algorithm [17]. To increase the reliability of the results we have used the mul- 
tiple histogram reweighting, due to Ferrenberg and Swendsen [18] along with the Lee and 
Kosterlitz's method [19] of finite size scaling for a first order phase transition. 



2 The definition of the thermodynamic quantities re- 
lated to the model 

The Monte-Carlo simulation was carried out on a square lattice of dimension L x L with 
the two-dimensional spins located at each site and interacting with the nearest neighbors 
via the Hamiltonian 

tf = £ 2 [l-(cos^i) p2 ] (2) 

(ij) 

The specific heat at a dimensionless temperature T is related to the energy fluctuation 

((H 2 ) - (H) 2 ) 

C v = ^ '- — )-LL 3 

NT 

where N is the number of spins. The conventional long range order parameter is given by 

(Pi) = (cos6) (4) 

where 9 is the angle that a spin makes with the preferred direction of orientation and the 
average is over the entire sample. The next higher rank order parameter is defined as 

(P 2 ) = 1 -(3cos 2 9-l) (5) 

The order parameter susceptibility is defined in terms of the fluctuations of the order pa- 
rameter (Pi) 

(6) 



The first rank pair correlation coefficient is denned as 



G 1 (r) = ((coaB li )) r 



(7) 



where % and j are two spins separated by a distance r. The second rank pair correlation 
coefficient is defined as 



3 The computational details 

In this section we briefly describe the Wolff cluster algorithm, the Ferrenberg-Swendsen 
multiple histogram reweighting technique and the Lee-Kosterlitz finite size scaling for first 
order phase transition. The Monte-Carlo (MC) simulations were performed on square lattices 
of size L 2 for L = 16, 32, 64, 96, 128, 160 and 192. We have used Wolff's cluster flip 
algorithm, the essential steps of which are as follows. 
(1) A random unit vector r is taken and a spin flip a x — > a' x is defined as 



G 2 (r) = (P 2 (cos^)) r 



(8) 



2(a x ,r)f 
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(2) Bonds (x, y) of the lattice are activated with a probability 
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where, 
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This process leads to the formation of a cluster on the lattice. 



(3) All spins in a cluster are now nipped according to a x — > a' x 

We have calculated the thermodynamic quantities using multiple histogram reweighting 
technique of Ferrenberg and Swendsen [18] , which is briefly described below. The partition 
function of the system is given by 

Z(^) = J2p(E)ex P [-(3E] (11) 

E 

where p(E) is the density of states, = 1/T (the Boltzmann constant has been set equal 
to unity) and E is the energy of the system. In the histogram reweighting method, energy 
histograms are generated at a number of temperatures with i = 1, 2, ....R and Ni(E) is the 
histogram count for the i th temperature. We denote by rij the total number of configurations 
generated in the i th simulation, i.e, rij = \^N{(E). According to references [18] and [20] . 

E 

the best estimate of the density of states, obtained after histogram reweighting, is given by 

R 

P (E) = - ir ^ (12) 

^gj'Z-'expl-^E] 

i=i 

where, gi — 1 + 2rj, Ti being the auto-correlation time for energy at the i th temperature. 
Substituting eqn.(13) in eqn.(12) gives us a self consistent equation for the partition function 
at any temperature (3: 

]T^iV.(£)exp[-/3£] 

e l^ g j n i Z i ex P\--PjE\ 
j 

One can also carry out the computation in terms of the probability instead of the parti- 
tion function. The unnormalized probability for an energy E in the k th simulation is given 
by 

Pk (E) = p(E)exp[-p k E] (14) 

i,e., 

E 

The free energy at the temperature (3 k is 

f k = -llnZ(^) (16) 

Pk 



i.e., 

exp(-&/ fc ) = = Z(f3 k ) (17) 

E 

In place of eqn.(14) we get a self consistent equation for the reweighted probability, 

Pk{E) = ■ (18) 

3 

The best estimate of the energy or any other observable Q is given by [20] 

(O) = 1 V 9i 1 Q ls eM-( 3E is\ (1Q) 

3 

where Qi S is the histogram count of the observable Q in the state s obtained during the i th 
simulation and E is is the total energy of such a state. The factors g k now correspond to the 
observable Q. Lee and Kosterlitz [Tj5] proposed a convenient method for the determination of 
the order of the phase transition which can be applied to systems having linear dimension less 
than the correlation length. For a temperature driven first order transition in a finite system 
of volume L d with periodic boundary condition one needs to compute the histogram count 
of the energy distribution denoted by N(E; f3, L). The p 2 = 50 model has a characteristic 
double peak structure for N(E; (3, L) in the neighborhood of the transition temperature. 
The two peaks of iV at E\(L) and E^{L) corresponding respectively to the ordered and 
disordered phases are separated by a minimum at E m (L). A free-energy-like quantity is 
defined as 

A(E; (3, L,N) = - In N(E; (3, L) (20) 

where M is the number of configurations generated. The quantity A(E; (3, L, A/") differs from 
the free energy F(E; (3, L) by a temperature and dependent additive quantity. A bulk 
free energy barrier can therefore be defined as 

AF(L) = A(E m , (3, L, A") - A{E l] (3, L, A") (21) 

It may be noted that at the transition temperature, A(E\\ p, L, J\f) = A(E 2 ; (3, L,J\f) and 
AF is independent of Af. For a continuous transition AF(L) is independent of L and for 
a temperature driven first order transition it is an increasing function of L, even when L is 



smaller than the correlation length, £, prevailing at the system at the transition temperature. 
If one is in a region where L is much greater than £, then AF obeys the scaling relation [19] 

AF ~ L d - X (22) 

Clearly, the temperature at which the double-well structure of A has two equally deep 
minima gives a precise estimation of the transition temperature. 

4 Results and discussion 

Square lattices of linear dimension L ranging from 16 to 192 were simulated and for each 
lattice simulations were performed at 9 to 13 temperatures in the neighborhood of the 
transition to record the histograms for energy. In Table [1] we have depicted for each lattice 
size and temperature the technical quantities of interest. These include the number of Wolff 
clusters generated (n c ), the percentage of average cluster size in units of lattice size for 
each temperature (c), the number of equivalent Monte Carlo sweeps (MCS) and the energy 
auto-correlation time r e . The number of configurations generated ranges from about 10 8 to 
10 9 . In FigJUwe have plotted the percentage of average cluster size in units of the lattice 
size (c) against temperature for L = 128. It is clear that the average cluster size for a given 
lattice, decreases with increase in temperature and there is a sharp fall at the transition. The 
maximum cluster size in units of the lattice size is about 84.4% (for L = 16) and is seen to 
decrease with increase in the system size. The auto-correlation time, which was calculated 
by the method proposed by Madras and Sokal [21], is seen to increase rapidly with the 
increase in lattice size and possesses a sharp maximum at the transition temperature. The 
logarithm of the peak value of the energy auto-correlation time has been plotted against L 
in Figj2] . We find empirically a scaling rule, lnr e ~ where = 3.05. The behavior of 
the order parameter correlation time r Q is found to be similar in nature as that of r e . 

The energy histograms obtained for L = 128 are shown in FigJHJ For this lattice, simula- 
tions were performed at 13 temperatures ranging from 1.0000 to 1.0175. This temperature 
range is rather small and were chosen to bracket the transition temperature. This diagram 
shows that there is an energy range where almost no sampling takes place for any temper- 
ature and there are dual peaked histograms at a number of temperatures where sampling 
takes place with one peak in the ordered phase and the other in the disordered phase. The 
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Figure 1: The average cluster size (c) in percents of the lattice size for L = 128 obtained 
using the Wolff algorithm. 
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Figure 2: Logarithmic plot of the peak value of the energy auto-correlation time against L. 

existence of these dual peaked histograms is a signature of a first order phase transition 
where two phases can coexist at a given temperature. 

The error in estimating the reweighted probability Pk{Q) from the raw histograms is 
given by 

Sp k (Q) = — ^Pk (Q) (23) 

n=l 

and this can be estimated directly from the histogram counts. The percentage error in 
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Figure 3: The histograms for E, the average energy per particle generated for the 128 x 121 
lattice for the p 2 = 50 model at the 13 temperatures indicated. 
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Figure 4: The average energy per particle E plotted against dimensionless temperature T 
for different lattice sizes. For clarity error bars are shown only for L = 16, 160 and 192. 

the reweighted probability for energy in the lattice L = 192 is about 0.74% where the raw 
histograms have peaks in the ordered phase. In the intermediate energy range where little 
sampling takes place for any choice of temperature the error is evidently large and this 
cannot be significantly reduced by any realistic effort. 

FigJH shows the temperature variation of the energy for a number of lattices, as is 
obtained by applying histogram reweighting technique. From the energy histograms, we 
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Figure 5: The free-energy barrier height AF plotted against lattice size L with the linear 
fit represented by straight line. 

have calculated the free energy like quantity A, defined as A(E; (3, L,Af) = — In N(E; (3, L) 
where N(E; (3, L) is the histogram count of the energy distribution. The free energy barrier 
AF(L) was evaluated and in Figj5]we have plotted AF against L where a good linear fit 
has been obtained. This is a direct verification of the scaling rule AF ~ L d ~ x of Lee and 
Kosterlitz [19] since the lattice dimensionality d = 2 in this model. We further note that 
the scaling relation is well obeyed down to L = 16 which happens to be of the order of the 
correlation length, £ for the system, as one can estimate from the relation AF (£) ~ 1 |19j . 

The specific heat C v was obtained from the energy fluctuation and Figj6] shows its tem- 
perature variation. It is evident that the peak height of C v grows rapidly at the transition. 
From FigJTJ where the maxima of C v are plotted it is clear that the standard scaling rules 
C v ~ L d for first order transition [22] are accurately obeyed in this model. 

We have also tested the the finite size scaling relation 

T c (L)-T c (oo)~L- d (24) 

which is valid for a first order phase transition [22]. T c (oo) represents the thermodynamic 
limit of the transition temperature T c . We have estimated the transition temperature in 
two ways - - T c " is the estimate of T c obtained from the peak position of the specific 
heat C v and T C F represents the transition temperature obtained from the fine tuning of the 
free energy vs energy curve to obtain two equally deep minima. In FigJS] the transition 
temperatures thus obtained have been plotted against L~ 2 . It is seen that the linear fits are 
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Figure 6: The specific heat C v plotted against temperature T for different lattice sizes. For 
clarity only the above lattice sizes are shown and the error bars have been indicated for two 
lattice size. 
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Figure 7: The peak heights of C v plotted against L 2 with the linear fit represented by the 
straight line. The error bars for most points are smaller than the dimensions of the symbols 
used for plotting. 

good within statistical errors and the thermodynamic limit of the transition temperature is 
1.00897 ± 6 x 10 -5 , within which the two linear fits are seen to converge. 

The pair correlation functions Gi(r) and G^r) were calculated for temperatures T = 
1.0081, 1.0085, 1.0092 and 1.0095 for L = 128 and are shown in FigJH and and \M The 
first two of these temperatures are less than the transition temperature for this lattice while 
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Figure 8: The transition temperature T c obtained from (a) specific heat peak position and 
(b) fine tuning of free energy curve plotted against L~ 2 along with the respective linear fits. 
The intercept on the Y-axis is 1.00897 ± 6 x 10~ 5 . 
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Figure 9: The plots of the pair correlation function G\(r) against r for the 128 x 128 lattice 
for the temperatures indicated. The curves are plotted for r ranging up to L/2. 

the other two temperatures are in the disordered phase. The curves have been fitted to a 
power law Gj(r) = air~ bi + fa for i — 1 and 2. It may be noted that the parameter / is 
the asymptotic value of the pair correlation function. We observe that while the first rank 
correlation function Gi(r) decays to zero at the two higher temperatures (f\ = 0), this is 
not the case for the higher rank correlation function G 2 (r) (f 2 ~ 0.22). In other words, 
while the lowest rank correlation among the spins vanishes just above transition, the next 
higher rank correlation continues to persist. 



0.8 
0.6 
0.4 

O 

0.2 


20 40 60 

r 

Figure 10: The plots of the pair correlation function G^fr) against r for the 128 x 128 lattice 
for the temperatures indicated. The curves are plotted for r ranging up to L/2. 

5 Conclusions 

The simulations in the two dimensional modified XY-model presented in this communication 
show that all the first order finite size scaling rules are obeyed. Computation has been 
performed in system size up to 192 x 192 which may be considered to be reasonably large for 
the purpose of arriving at a conclusion regarding the behavior of the model. We are inclined 
to conclude that the model exhibits a first order phase transition. This is in agreement 
with the views of some of the earlier investigators including Domany et.al [1] and van Enter 
and Shlosman [IB] . The existence of a quasi-long-range-order- disorder transition observed 
in the 2-D XY-model is known to be due to vortex-antivortex unbinding (KT transition). 
In absence of the role played by the vortices, one would not observe any order- disorder 
transition in the XY-model in accordance with the Mermin- Wagner theorem. In the class 
of models we have investigated the role played by the vortices changes qualitatively with 
change in p 2 (which increases the non-linearity of the potential well) as has been seen in 
the early work of Himbergen [6] . Also we have seen that the number of vortex pairs grows 
rapidly with the increase in p 2 [23J . Qualitatively, one may therefore think that the modified 
XY-model for large values of p 2 , behaves like a dense defect system and gives rise to a first 
order phase transition as has been predicted by Minnahagen [HI QUI ED] • 

A similar change in the nature of phase transition has been observed to occur in a two- 
dimensional Lebwohl-Lasher model and a modified version of it [21]. The potentials in these 
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two models are — P 2 ( C0S 6 I ) and — P 4 (cos6') respectively, the latter having a greater amount 
of non-linearity. Although both models are in the same universality class, it was observed 
that while the — P2(cos8) potential leads to a continuous transition, the modified model 
with — P 4 (cos#) potential exhibits a strong first order phase transition. It has also been 
noticed that the suppression of the defects in these models leads to a total disappearance 
of the phase transitions [25] . 

We mention another point before ending this section. This is the performance of Wolff 
cluster algorithm which turned out to be very convenient to simulate the model. Conven- 
tional algorithms, as we have seen, does not work well in this model. Our earlier attempt 
[T6] using the recently developed Wang-Landau (WL) algorithm [26] which directly deter- 
mines the density of states of a system is also not a good choice for simulating this model. 
The main problem while using the WL algorithm is that configurations near the minimum 
energy take a very long time to be sampled during the random walk and it becomes imprac- 
tical to simulate continuous models of even moderate size because of the huge CPU time 
that becomes necessary. Among other things, a great virtue of the Wolff algorithm is that 
it does not contain any adjustable parameter even while simulating a continuous model. 

Besides using the Wolff algorithm for the simulation we have used the Ferrenberg- 
Swendsen multiple histogram reweighting technique and the finite size scaling rules of Lee 
and Kosterlitz. We conclude by noting that a combination of these computational tools 
till now provides a very efficient and accurate method of analyzing results obtained in an 
unknown system. 
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Table 1: T is the dimcnsionless temperature, L is the lattice size, n c is the number of Wolff clusters , 
which varies from 1.1 x 10 8 to 10 9 with lattice sizes, (c) is the average cluster size as percent of the lattice 
size, MCS is the number of equivalent Monte Carlo sweeps in units of 10 8 and r e is the auto-correlation 
time for energy (in units of Wolff clusters). 
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5.60 


1.81 


0.09 


0.049 


0.031 


0.023 










Te 


239 


354 


788 


3873 


1851577 


54589 


44409 


37398 


39685 
























L = 128 














T 


1.0000 


1.0025 


1.0050 


1.0075 


1.0080 


1.0085 


1.0090 


1.0092 


1.0095 


1.0100 


1.0125 


1.0150 


1.0175 


n c 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


10 9 


(c) 


69 


67 


65 


61 


59 


58 


41 


9 


1 


0.76 


0.36 


0.25 


0.19 


MCS 


6.96 


6.78 


6.55 


6.14 


5.99 


5.82 


4.17 


0.92 


0.11 


0.076 


0.036 


0.025 


0.019 


T e 


193 


242 


344 


780 


1111 


1906 


838017 


814172 


33708 


27202 


24802 


21925 


20478 












L = 96 


















T 


1.0000 


1.0025 


1.0050 


1.0075 


1.0080 


1.0085 


1.0090 


1.0093 


1.0100 










n c 


9 x 10 s 


9 x 10 8 


9 x 10 8 


9 x 10 8 


9 x 10 8 


9 x 10 8 


9 x 10 8 


9 x 10 8 


9 x 10 8 










(c) 


70 


67 


65 


62 


60 


58 


47 


20 


1 










MCS 


6.30 


6.14 


5.93 


5.58 


5.45 


5.23 


4.25 


1.86 


1.33 










T e 


185 


233 


335 


658 


3720 


49156 


291844 


139640 


20637 






















L = 64 
















T 


1.0000 


1.0025 


1.0050 


1.0075 


1.0081 


1.0087 


1.0093 


1.0100 


1.0125 


1.0150 


1.0175 


1.0200 






8 x 10 s 


8 x 10 8 


8 x 10 8 


8 x 10 8 


8 x 10 8 


8 x 10 8 


8 x 10 8 


8 x 10 8 


8 x 10 8 


8 x 10 s 


8 x 10 8 


8 x 10 s 




(c) 


70 


69 


66 


63 


61 


56 


43 


15 


1 


1 


0.78 


0.63 




MCS 


5.66 


5.52 


5.34 


5.04 


4.92 


4.55 


3.46 


1.27 


0.11 


0.08 


0.06 


0.05 




T e 


181 


229 


319 


977 


2451 


37018 


61038 


54839 


6548 


5755 


5714 


5519 
















L = 32 
















T 


0.9900 


1.0000 


1.0050 


1.0075 


1.0100 


1.0112 


1.0125 


1.0137 


1.0150 


1.0200 


1.0300 


1.0400 




n c 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 8 


1.7 x 10 s 


1.7 x 10 8 


1.7 x 10 s 




(c) 


76 


72 


68 


65 


52 


39 


23 


11 


7 


2 


1 


1 




MCS 


1.29 


1.23 


1.16 


1.10 


0.88 


0.67 


0.39 


0.19 


0.11 


0.04 


0.02 


0.019 




T e 


110 


170 


561 


2149 


7488 


8420 


7941 


6152 


4070 


1556 


1355 


1290 


















L = 16 














T 


0.9500 


0.9750 


1.0000 


1.0062 


1.0125 


1.0188 


1.0219 


1.0250 


1.0312 


1.0375 


1.0500 


1.0800 


1.1000 


n c 


1.1 x 10 8 


1.1 x 10" 


1.1 x 10 8 


1.1 x 10 8 


1.1 x 10 8 


1.1 x 10 8 


1.1 x 10 8 


1.1 x 10 8 


1.1 x 10 8 


1.1 x 10 s 


1.1 x 10 8 


1.1 x 10 s 


1.1 x 10 8 


(c) 


84 


80 


74 


66 


57 


38 


25 


17 


8 


5 


3 


2 


1 


MCS 


0.93 


0.889 


0.816 


0.734 


0.631 


0.425 


0.275 


0.197 


0.094 


0.060 


0.042 


0.025 


0.021 


Te 


73 


82 


269 


657 


1226 


1593 


1560 


1259 


736 


507 


332 


319 


281 



